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1. Introduction 


ABSTRACT 





Recently, a coupled echo canceller was proposed that uses two short adaptive filters for 
sparse echo cancellation. The first filter operates in the partial Haar domain and is used 
to locate the channel’s dispersive region; the second filter is then centered around this 
location to cancel the echo in the time domain. In this paper, we propose feasible 
solutions to improve the performance of this partial Haar dual adaptive filter (PHDAF) in 
practical applications. These include: (1) alleviating the dependence of the PHDAFs 
performance on the echo-path impulse response’s bulk delay; (2) improving the 
tracking performance of the PHDAF in response to abrupt changes in the echo path; and 
(3) extending the original PHDAF structure to support the cancellation of multiple 
echoes. The proposed algorithmic solutions exploit the Haar transform’s polyphase 
representation and make use of a novel peak tendency estimator (PTE) based on 
Dezert-Smarandache theory (DSmT). The improved PHDAF is evaluated in terms of its 
mean-square error (MSE) curves and its mean time to properly locate a dispersive region 
for different SNRs. Results show that enhanced performance can be obtained using the 
proposed solutions at a minimal increase in computational cost. 

© 2008 Elsevier B.V. All rights reserved. 





severely impede a conversation. The coding and signal 
processing functions of digital technologies may introduce 


Line, or network echo is commonplace in today’s 
expanding communications infrastructure. Unlike other 
types of echo (e.g. acoustic), line echo is sparse: the echo- 
path impulse response consists of an initial zero or bulk 
delay region, corresponding to the signal round-trip, 
followed by a non-zero or dispersive region, corresponding 
to the echo arrival. Line echo is usually caused by an 
HOPET mismatch, as occurring in the hybrid circuits 
used for 4wire conversion [1]. In voice communications, 
when the bulk delay between callers exceeds 25 ms or so, 
the reflected signal is perceived as a distinct echo that can 
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delays in excess of 100ms; while for long distance calls 
routed via satellites, the propagation delay may reach 
several 100 ms [2]. Recent advancements, such as Voice 
over Internet Protocol (VoIP) telephony and xDSL tech- 
nologies for broadband data transmission, highlight the 
need to develop better echo cancellers for sparse line 
echo. 

In one of the earliest works on sparse echo cancellation 
[3], the input and desired signals are bandpass-filtered, 
decimated, and used by a short adaptive filter to estimate 
the bulk delay. A second short filter operating at the 
original sampling rate is centered around the dispersive 
region to cancel the echo. This way, only two short 
adaptive filters are required (compared to one long filter), 
thus reducing overall system complexity. However, the 
use of bandpass filtering in this approach has important 
drawbacks. Firstly, in speech applications, it can remove 
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important frequency components in the signals, prevent- 
ing proper convergence of the subsampled adaptive filter. 
Secondly, bandpass filtering in effect smears the peak of 
the unknown echo-path impulse response, making it 
more difficult to correctly locate the dispersive region. 

Recent literature is rich in adaptive filtering algorithms 
that exploit the sparse characteristics of line echo [4]. 
Most of these algorithms are based on finding ways 
to determine which filter coefficients are associated 
with the echo, and then adapting only these coefficients. 
An adaptive multiple echo (ME) canceller is proposed in 
[5], which uses a full-length primary adaptive filter in 
parallel with a group of short secondary adaptive filters. 
In [6], a two-stage adaptation process is proposed in 
which the first stage estimates the bulk delay while the 
second stage adapts filter coefficients using a constrained 
tap-selection approach [7]. A well-known class of sparse 
echo cancellers are based on the proportionate normal- 
ized least mean squares (PNLMS) algorithm [8] and its 
variants [9-11]. These algorithms allocate individual step- 
size gains in proportion to the magnitude of each filter 
coefficient. 

Other recent attempts at improving sparse echo 
cancellation rely on applying orthogonal wavelet trans- 
forms to the input data. It is shown in [12] that the 
number of adaptive coefficients needed to cancel the echo 
can be reduced significantly by applying a Haar transform. 
In [13], the authors propose using a subset of Haar 
wavelets to detect the significant channel coefficients. By 
exploiting the hierarchical structure of the dyadic wavelet 
expansion, these significant coefficients are used to 
activate wavelets in the remaining Haar subsets that 
share the same non-zero time-support. The Haar trans- 
form is particularly attractive for sparse echo identi- 
fication: it facilitates the location of the dispersive 
region through proper selection of its scale/transla- 
tion parameters, and it is easily amenable to a digital 
implementation. 

Recently, Bershad and Bist [14] have proposed a 
solution to the sparse echo cancellation problem that 
combines favourable attributes of [3,13] in a coupled 
configuration consisting of two short adaptive filters. In 
this approach, referred to here as the partial Haar dual 
adaptive filter (PHDAF), the first filter operates on a subset 
of input Haar coefficients, and is used by a peak delay 
estimator to locate the echo-path’s dispersive region. The 
second filter is centered around this location to cancel the 
echo in the time domain. In cases where the bulk delay is 
large, the PHDAF provides a significant reduction in 
computational and memory requirements. In addition, 
by reducing the number of filter taps to an amount 
necessary to model the dispersive region, the convergence 
speed of the echo canceller is increased. An improved 
theoretical model of the LMS algorithm in [14] for low 
rank systems was recently developed in [15] to better 
predict the behaviour of first and second moment 
Statistics of the partial Haar adaptive filter. In [16], a 
partial block wavelet transform is proposed to increase 
efficiency and improve peak detection of the system in 
[14]. It is shown that the block transform reduces to using 
Daubechies’ biorthogonal 2.2 spline wavelet, which is 


claimed to have better properties for estimating a peak’s 
location. 

In this paper we identify, and propose feasible 
solutions to three inherent limitations of the PHDAF for 
Sparse echo cancellation in [14]: (1) dependence of the 
PHDAFs performance on the echo-path impulse response’s 
bulk delay; (2) degraded tracking performance of the 
PHDAF in response to abrupt changes in the echo-path 
impulse response; and (3) limitation of the original 
PHDAF to a single dispersive region. The proposed 
algorithmic solutions exploit the polyphase representa- 
tion of the Haar transform and make use of a novel peak 
tendency estimator (PTE) based on Dezert-Smarandache 
theory (DSmT) and fuzzy inference [17,18]. The improved 
PHDAF is evaluated in terms of its mean-square error 
(MSE) curves as well as its mean time to properly locate a 
dispersive region under different SNRs. In experiments 
using normalized least mean squares (NLMS) for the filter 
coefficient adaptation, the proposed amendments to the 
original PHDAF are shown to yield significant performance 
gains at a minimal increase in computational cost. 

This paper is organized as follows: The structure and 
main equations of the PHDAF in [14] are reviewed in 
Section 2, along with a discussion of its main limitations. 
The solutions that we propose to overcome the latter are 
developed in Section 3. A series of supporting computer 
experiments is presented in Section 4. Finally, Section 5 
concludes the work. 


2. Background and problem formulation 
2.1. The partial Haar dual adaptive filter 


2.1.1. Partial Haar transform 

Let N=2/, where J is a positive integer. The 
N-dimensional discrete-time Haar wavelet transform can 
be represented by an N x N orthogonal matrix H with 
entries hm(n), where the row index m identifies a basis 
vector and the column index n represents discrete-time 
[19]. The elements of the first row are equal to 
ho(n) = 1/VN, while the remaining rows are obtained by 
scaling and shifting a discrete-time wavelet filter y(n) 
defined as 


+1, O<n<N/2, 
W(n)=<« —-1, N/2<n<N, (1) 
0 otherwise. 


Specifically, let m = 2 + k, where j e {0,...,J—1} is the 
scale index and k e {0,...,2/ — 1} is the translation index. 
We have 


hm(n) = ajw(2!n — kN), (2) 


where the normalization factor &j = \/ 2! /N. The number 
of rows with scale index j is equal to X and the 
corresponding basis vectors have non-overlapping time- 
support of length N /2) . A partial Haar wavelet transform 
consists of using only a subset of size q = 2! of Haar basis 
vectors corresponding to scale index j, to transform a 
given input data vector. The corresponding transform 
matrix, denoted H4, is thus a qx N submatrix of H, 
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obtained by extracting its rows h,,(n) for m = q + k and 
k=0,...,q—1. As j increases, the time resolution of the 
transform improves at the cost of increasing the dimen- 
sion of the output vector, which translates into a longer 
Haar-domain adaptive filter in the PHDAF (see below). 
Unlike the complete Haar transform matrix, the partial 
Haar matrix Hg (q<N) is only row-wise orthonormal, i.e. 
HH, = 1, 4 H; H4. 


2.1.2. PHDAF structure and operation 

The structure of the PHDAF proposed in [14] is shown 
in Fig. 1, where discrete-time signals u(n) and d(n), 
respectively, denote the line input and the echo-contami- 
nated reference used to drive the adaptive process. 

The upper branch, consists of a partial Haar transform 
matrix H4 of size q x N, followed by a length-q (<N) partial 
Haar adaptive filter. The value of N is set to match the 
maximum length of the unknown echo-path impulse 
response. The selection of q = Z involves a trade-off between 
the desired level of time resolution and complexity reduction. 
At time n, a new sample u(n) is shifted into the input data 
vector u(n) = [u(n), u(n — 1),..., u(N — 1)]! of length N and a 
new transformed input vector of length q is calculated, i.e. 


z(n) = Hgu(n). (3) 


Vector z(n) is used as input into a so-called partial Haar 
adaptive filter with coefficient vector at time n given by v(n). 
In this work, the partial Haar adaptive filter v(n) is updated 
by means of the NLMS algorithm, i.e. 


v(n + 1) = v(m) + pullz(n)||-*en(1)z(n), (4) 


ey(n) = d(n) — v(n)'z(n), (5) 


where u denotes the step size and e,(n) is the error signal. 
Under appropriate conditions, v(n) converges to a scaled- 
down version of the true echo-path impulse response. The 
peak delay estimator in Fig. 1 tracks the location of the 
dispersive region by locating the peak coefficient magnitude 
of v(n). 


In the lower branch of the PHDAF, the estimated peak 
location is used to offset a short time-domain adaptive 
filter w(n) of length L so as to center it around the 
dispersive region. The value of L is set to match the longest 
expected dispersive region in the echo path. Once the 
peak location of v(n) is properly estimated, the short time- 
domain adaptive filter can converge to the dispersive 
region and thus cancel the echo component of the 
reference signal d(n). In this work, the NLMS algorithm 
is also used to update w(n). In practice, the residual signal 
after echo removal, i.e. e(n) = d(n) — w(n)'u(n), is used for 
data transmission over the network. 


2.1.3. Properties of the PHDAF 
To simplify the discussion, let us assume that 


d(n) = w,u(n) + v(n), (6) 


where Wo represents the true echo-path impulse response 
and wn) is an additive measurement noise. Also assume 
that the sequences u(n) and v(n) are stationary white, 
uncorrelated, with variances o2 and o2, respectively. Then, 
the optimum Wiener solution for the partial Haar filter 
can be obtained as 


Vo = H,Wo, (7) 


which is the partial Haar transform (i.e. a scaled-down 
version) of the true time-domain echo-path impulse 
response Wo. Thus, when Wo features a single, localized 
dispersive region and the scale index j=log,q of the 
partial Haar transform provides adequate temporal reso- 
lution, the location of the dispersive region can be 
recovered from the peak coefficient in Vo. It is also worth 
mentioning that Wo cannot be recovered from vo unless it 
lies in the row-space of H4. Under proper operating 
conditions, the partial Haar adaptive filter will converge 
to the Wiener solution in (7). The steady-state MSE then 
takes the form J =Jmin + -@W where Mæ is the misadjust- 
ment noise power of the adaptation algorithm and Jmin is 
the minimum MSE corresponding to the Wiener solution 





Fig. 1. PHDAF structure. 
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(7). The latter is given by 


Jmin = 07 + 02(||Woll5 — IIVoll5)- (8) 


When q = N (complete Haar transform), we have ||Wo|l2 = 
[Vo llz and the second term on the right-hand side of (8) is 
equal to zero. When q<N, this term is positive, reflecting 
the rank-deficiency of the transform matrix Hg. This 
increase in minimum MSE is acceptable if it does not 
impact the ability of the partial Haar adaptive filter to 
locate a dispersive region. 

These considerations motivate the use of the con- 
verged vector v(n) as input to the peak delay estimator for 
the purpose of locating the dispersive region. It is shown 
in [14] that the greater the steady-state peak magnitude of 
v(n), the faster the short time-domain filter can be 
centered, and therefore, the overall PHDAF can adapt 
more quickly. 


2.2. Limitations of the PHDAF 


The simulation results in [14] show that the PHDAF 
may provide a drastic increase of the echo canceller’s 
convergence speed in general, while keeping computa- 
tional complexity low. However, our own experiments 
show that in certain cases, the PHDAF can require an 
extremely large number of input samples to converge, 
underlining inherent shortcomings of the original system. 
Below, we take a closer look at these drawbacks. 


2.2.1. Sensitivity to bulk delay 

The lack of shift-invariance of the Haar wavelet 
transform can greatly affect the time required by the 
peak delay estimator to properly locate a dispersive 
region. A wavelet transform whose basis vectors have a 
time-support of length M = N/q is periodically shift- 
invariant with a period equal to M. As a result, depending 
on the bulk delay of the true echo-path impulse response, 
there exist M distinct partial Haar transformed Wiener 
solutions. Furthermore, the respective peak magnitudes of 
each of these solutions can differ, and consequently, the 
amount of time it takes the peak delay estimator to 
correctly locate a dispersive region. 

To illustrate this point, Fig. 2 shows the transformed 
impulse responses of differently shifted versions of a 
common echo-path impulse response. In this example, the 
true echo-path vector of length N=1024 is set to 
Wo = [m(n — 600 — 1)]', where m(n) is one of the ITU-T 
G.168 hybrid impulse responses in [20] and shift index 
l e {1,5,6,7} (see Section 4 for additional details). For each 
value of l, a transformed impulse response vector of size 
q = 128 is computed by applying the partial Haar matrix 
H; as in (7). One can clearly observe that the second case 
(l= 5) displays a diminished peak. 

Due to measurement noise and the rank-deficiency 
effect of the partial Haar transform on the minimum 
MSE in (8), the peak delay estimator takes a longer time 
to locate such small peaks. To support this claim, consider 
the family of experimental cumulative distribution 
function (CDF) curves in Fig. 3 corresponding to the 
time! to correctly locate the dispersive region for these 





0 20 40 60 80 100 120 


samples 


Fig. 2. Set of transformed impulse responses corresponding to m(n — 
600 — D for shift index l € {1,5,6,7}. 





10° 104 


No. of Iterations 


Fig. 3. Probability (P.) of correctly estimating peak delays within a given 
number of time iterations (each curve is labelled by its peak magnitude 
in the Haar domain and the corresponding value of 1). 


four different bulk delays. It is apparent from this figure 
that it takes the peak delay estimator a larger number of 
input samples to correctly locate a smaller peak. The curve 
corresponding to a peak magnitude of 0.0166 shows that 
the peak delay estimator can barely locate the peak after 
107 input samples. 





1 The peak delay estimator looks for the filter coefficient with the 
maximum instantaneous absolute value, and returns its location A(n). If 
|A(n) — Ag|<y for n: consecutive input samples starting with n = nqa, 
where Ao is the true peak delay, then ng is set as the time to correctly 
estimate the peak location. Integer 7 is chosen to avoid small estimate 
jitters around Ao. As pointed out by Ribas et al. [16], a similar solution to 
the jitter problem was used in [14]. In Fig. 3, we set 7 = 10, and n: = 950. 
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2.2.2. Degraded tracking performance 

In [14], the merits of using the coupled PHDAF 
configuration are demonstrated by calculating the echo 
canceller’s mean time to correctly locate a dispersive 
region. It is shown on average that this time does not 
exceed more than 200 samples at 8kHz when employing 
a partial Haar NLMS algorithm. These results, however, 
only represent the case when the partial Haar adaptive 
filter has been initialized to zero, i.e. v(0) = 0, and the 
echo path is stationary. The case where the echo-path 
impulse response changes during the course of the echo 
canceller’s steady-state operation is not considered. Our 
experiments show that the performance of the PHDAF can 
degrade significantly in the presence of sudden changes in 
the true echo path. Indeed, unlike the case where v(0) = 0, 
and a peak only has to compete with neighbouring 
coefficient noise, after a sudden change in bulk delay, 
the new peak also has to compete with the residual of an 
older peak. In the case of the LMS adaptive algorithm for 
example, this older peak will undergo a slow exponential 
decay which may prevent the detection of the new peak 
for a large number of time iterations. 


2.2.3. ME cancellation requirements 

The peak delay estimator analyzed in [14] assumes the 
presence of a single dispersive region. This greatly 
simplifies the estimation process because the adaptive 
filter coefficient associated with a global maximum 
directly indicates the location of this region. In a multiple 
peak delay estimation scenario, the task does not only 
involve finding local maxima from the partial Haar 
adaptive filter, but also requires the classification of these 
maxima as peaks. A major difficulty arises from the fact 
that knowledge about the number of dispersive regions in 
a channel may not be available a priori. Clearly, this 
number must remain small if any advantages are to be 
gained from the use of a tailored ME canceller structure. 

Another difficulty relates to the shift-variant property 
of the Haar wavelets. It was seen that, depending on the 
bulk delay, the partial Haar Wiener solutions can display 
peaks of varying magnitude, and this in turn determines 
how fast a dispersive region is located. The same is true 
for MEs, except now their corresponding peaks with 
largest magnitude need not correspond to the same shift 
within a period of M. Finally, an algorithm that cancels 
MEs will require a dynamic mechanism to allocate 
resources (e.g. short time-domain filters) to dispersive 
regions that have been correctly located. 


3. Improving the coupled echo canceller 
3.1. Redundant partial Haar transform 


From the above discussion, it appears that some of the 
difficulties related to the rank-deficiency of the partial 
Haar transform, especially the sensitivity of the peak 
magnitude in the Wiener solution Vo (7) to the bulk delay 
in the true echo-path impulse response Wo, could be 
overcome by driving the partial Haar adaptive filter 
towards a properly shifted version of Wo. In effect, this 


can be achieved by modifying the upper Haar-domain 
adaptive filter in Fig. 1 in two different ways (without 
directly affecting the operation of the lower time-domain 
adaptive filter): 


(1) Replacing the reference signal d(n) in (5) by d(n+ D, 
where |! is the desired integer shift. 

(2) Replacing the transformed input vector z(n) by z(n — I) 
in (4)-(5). 


Under the modelling assumptions presented at the 
beginning of Section 2.1.3, both approaches lead to the 
Same modified optimum Wiener solution, i.e. 


Vo = H,S/Wo, (9) 


where S is an N x N shift matrix with (i,j)th entry equal to 
1 for j = i + l (i.e. Ith diagonal) and O otherwise. 

While both approaches have the equivalent effect of 
decreasing the effective bulk delay of Wo by | samples, in 
this work, we favour the use of the second approach. 
Indeed, as explained below, all M = N/q shifted input 
vectors z(n—1), where le {0,1,...,M—1}, can be ob- 
tained naturally as the polyphase components of a so- 
called redundant partial Haar transform (RPHT), which 
trades arithmetic operations for memory. In addition, 
shifting d(n) cannot help when dealing with MEs since a 
ME path channel can be written as a sum of single-echo- 
path channels, each with its own bulk delay, and whose 
contribution to d(n) is unknown a priori. 

A consequence of the wavelet transform’s lack of shift- 
invariance is the requirement of calculating the partial 
Haar transform of the length-N input data block u(n) in (3) 
for every new input sample. The computational load for 
this operation is on the order of N arithmetic operations 
per iteration,” which is expensive. The length-N vector of 
RPHT coefficients at time n, denoted s(n), is defined by 


s(n) = [s(n), s(n — 1),... s(n — N + DF, (10) 
where 
s(n) = h'u' (n), (11) 


h = &; [W(0), W(q),... W(q(M — 1))]' is a vector of length M 
containing the non-zero portion of any row of Hg, and 
u(n) = [u(n), u(n — 1),...,unm—M-+1)]'. Because this 
RPHT has a simple translational relationship with the 
input, it only requires the calculation of a single partial 
Haar coefficient s(n) every iteration with a small cost on 
the order of M arithmetic operations. 

To extract the q coefficients corresponding to the 
standard partial Haar transform z(n) = H,u(n) from s(n), 
we simply observe that 


z(n) = [s(n), s(n — M),...,s(n — (q — 1)M)]!. (12) 


In effect, vector z(n) corresponds to the first (i.e. l = 0) 
polyphase component of the complete vector s(n), while 
the remaining polyphase components represent the 





2 That is, N — q additions and q multiplications if the normalization 
factors «; in (2) are included. 
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Fig. 4. Block diagram: peak tendency estimator. 


Decreasing Fuzzy Rule Base 
Decreasing Peak 
transformed input vectors time delayed by !=1,...,M — 


1 sample(s), i.e. z(n — 1) down to z(n — M + 1). 

Unlike the recently proposed partial block wavelet 
transform-based algorithm in [16] which calculates a 
transformed input coefficient every M samples, the 
proposed approach here can produce a new transformed 
input sample s(n) every iteration. This makes it possible to 
drive the partial Haar adaptive filter with a specific 
polyphase component of s(n), which in turn results in its 
convergence to one of M different steady-state solutions 
(see Fig. 2, for example), each with its unique peak 
magnitude. 

In what follows, the term context will be used to denote 
the polyphase component of s(n), corresponding to shift 
index | € {0,...,M — 1}, that is being used as input to the 
partial Haar adaptive filter. Storing all N samples of s(n) 
makes M contexts and their resulting steady-state solu- 
tions readily available. At any given time, the partial Haar 
adaptive filter can be driven to one? of these M solutions, 
increasing the system’s flexibility by making it possible to 
avoid suboptimal cases (such as illustrated by the right- 
most curve in Fig. 3). 


3.2. Peak tendency estimation 


Without further a priori knowledge (e.g. a pre-deter- 
mined threshold), the magnitude of a detected peak alone 
is not a sufficient indicator of how well the bulk delay 
estimation process is progressing for a given context. 
A method is needed that includes a relative measure of 
the peak quality with regard to the detection process, in 
addition to a means for tracking this measure’s behaviour 
over time. To categorize the performance of a context, i.e. 
if it is suitable or not for quick and correct peak delay 
estimation, a PTE is therefore presented below. This 
estimator is based on the work of [21], where an original 
approach for target behaviour tendency estimation 
(i.e. receding or approaching) is developed using the 
DSmT of plausible and paradoxical reasoning [18]. 





3 An obvious, but much more expensive solution is to use M partial 
Haar adaptive filters in parallel, each operating with a different context 
and choosing the peak delay to correspond to the maximum peak among 
the global maxima for each filter. 


To mitigate the local effects of coefficient noise on the 
peak location of the partial Haar adaptive filter v(n), a 
peak discernibility measure (PDM) is proposed here which 
is calculated as follows: partition vector v(n) into three 
contiguous groups of similar size and find the maximum 
peak magnitude for each of the three groups. Let C,j,(n) 
and Cmax(n), respectively, denote the minimum and 
maximum among the three maxima so obtained. The 
PDM at time n is defined as* 


PDM(n) zls Cmin(n)/Cmax(N). (13) 


When Cmax(n) is large compared to Cmin(n), then the PDM 
approaches 1 (a large localized peak). However, if the two 
are comparable (which is usually the case for coefficient 
noise), then the PDM approaches 0. In Appendix A, we 
show that the PDM (13) is indeed related to the 
probability of correctly locating the peak of v(n), and thus 
provides an indication of the quality of the peak detection 
process. 

The proposed PTE takes PDM(n) as input at time n and 
uses a combination of DSmT and fuzzy inference similar to 
[21]. Unlike traditional logic, fuzzy inferences can be made 
even when rules are only partially satisfied. This is in stark 
contrast to crisp logic, where a rule is satisfied only if a 
premise of interest matches an implication’s antecedent 
exactly. This is a powerful extension for dealing with fuzzy 
sets as linguistic descriptors that are difficult to clearly 
define, and possess a degree of uncertainty (interested 
readers may consult [17,22-25] for a good introduction to 
fuzzy systems). This inference mechanism makes fuzzy 
systems well-suited for complex or nonlinear control 
applications, such as step-size control in adaptive filtering 
algorithms [26-28]. However, unlike such approaches, the 
method described here is solely used to track and 
categorize a peak’s discernibility as increasing or decreas- 
ing (where this information is used in the algorithm 
described in Section 3.3). 

The adopted PTE system is represented in block 
diagram form in Fig. 4. It consists of two state prediction 
filters running in parallel on the sequence of fuzzified 
inputs PDM(n), with each filter using a different model of 
possible PDM behaviour. The top model corresponds to a 
PDM that increases over time and eventually results in the 


4 Set PDM(n) = 0, if Cmax(n) = O. 
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Fig. 5. Fuzzy interface including small and large membership functions. 


proper detection of a dispersive region. The bottom model 
corresponds to a decreasing PDM which can be used to 
detect coefficient noise. The basic components of the PTE 
system are briefly described below (see [18,21,29] for 
additional details on DSmT and its applications to 
behaviour tendency estimation). 

Two fuzzy sets, labelled S and L and defined in terms of 
the membership functions shown in Fig. 5, are used to 
characterize a small and large discernibility, respectively. 
At time n, a new (row) vector of fuzzy PDM values is 
obtained as vn = [Vn(S), Vn(S U L), Vn (L), Vn(S A D] where 


Vn(L) = u, (PDM(n)), (14) 
Vn(S) = Us(PDM(n)), (15) 
Vn(S UL) = 1 — Va(S) — vn (L) (16) 


and va(SN L) = 0. The terms v,(-) denote the basic belief 
masses (bbm) used within the DSmT framework. 

The peak tendency at time n is internally characterized 
by a pair of fuzzy PDM state vectors corresponding to the 
two models, i.e. 4®° and pec, where (M = inc or dec) 


Hf = (ut (S), u (SU D, wf D, u (SOD), (17) 


with u“ (>) representing the corresponding bbm values. 
The state vector of each model is updated recursively in 
two steps: (1) the previous state vector u} is processed 
by a fuzzy inference engine to obtain a prediction of the 
current state vector, denoted tiai (2)the predicted state 
vector Mfp is merged with the current observation 
vector vn via the DSmT rule of combination to obtain the 
current state vector pw”. 

The respective rule-bases for each model in Table 1 
consist of if-then statements with antecedents A; and 
consequents C; (e.g. for rule 1 of the increasing PDM 
model, A; = C4 = S). These rules are used together with 
the membership functions in Fig. 5 to construct a 4 x 4 
fuzzy graph matrix G“ for each of the two models, as 
given in Table 2, which quantifies the degree of associa- 
tion between two fuzzy sets for a specific model governed 
by a corresponding rule-base. The graph entries for a 
transition A —> C are computed with the maximum value 


Table 1 
Fuzzy rule bases. 


Rule no. Increasing PDM 

1 If PDM(n — 1) = S then PDM(n) = S 

2 If PDM(n — 1) = S then PDM(n) = L 

3 If PDM(n — 1) = L then PDM(n) = L 
Decreasing PDM 

1 If PDM(n — 1) = L then PDM(n) = L 

2 If PDM(n — 1) = L then PDM(n) = S 

3 If PDM(n — 1) = S then PDM(n) = S 

Table 2 


Fuzzy graphs corresponding to two models of peak discernibility. 


m= il c= M S SUL L SOL 


(a) Increasing discernibility, G”: 


S 1 0 1 0 
SW 0 0 0 0 
L 0.2 0 1 0 
Sime 0 0 0 0 
(b) Decreasing discernibility, G% 

S 1 0 0.2 0 
Sie 0 0 0 0 
L 1 0 1 0 
SAL 0 0 0 0 


of the fuzzy intersection function between A and all A; 
such that C; = C [22]: 


Sasc = Max MIN(L,4(X), Ha (X)), (18) 
j.x¢[0,1] J 


where A and C take values in the DSmT set of composite 
propositions 2 = {S,L,S N L,S UL}. As long as the member- 
ship functions u4 and Ha, overlap, the degree of associa- 
tion is non-zero. The predicted state vector for each model 
is obtained as 


ae = we 0G”, (19) 


where o denotes the matrix form of Zadeh’s max-min rule 
of composition [17]. Each vector Ba is then normalized 
so that its entries sum up to 1. 

At time n, the new input vector vn is separately 
combined with each of the predicted state vectors mfn 
using the DSmT rule of combination 


O= X hinm A B, CE2 (20) 
A,BE2, 
ANB=C 
to obtain the updated fuzzy state vectors pw”. 
Finally, each model’s corresponding set of updated 
bbms (u) are converted into the Pignistic probabilities 
[18] 


P {S} = w (S) + 0.54” (SU D +0.5u“(S NL), (21) 


PL} = w“ (L) + 0.5u“(SUL) + 0.5u“(S OL). (22) 
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The correct behavioural model corresponds to the model 
with the smallest Pignistic entropy, as given by 


He =— X. P“{A}In(P“ {A)). (23) 
Ac{S,L} 


3.3. Suboptimal-context escape (SCE) 


The SCE approach that we propose here monitors the 
correct peak model from the PTE and decides if the 
current context should be changed (i.e. if a different 
polyphase component of the input RPHT should be used) 
in the event that a peak’s discernibility is poor. Deciding 
when to make a context switch requires some care. While 
the partial Haar adaptive filter needs time to reach steady 
State, it is usually the case that the time required to 
correctly locate a peak is far less than the convergence 
time in a proper operating context. It is thus inefficient to 
wait until the Haar-domain adaptive filter has converged 
to decide whether to remain in the current operating 
context or to attempt a new context. 

The proposed solution consists of a schedule of M = 
N/q  non-decreasing trial periods t= {T,...,Ty} 
(T1 <T2<--- <Ty) to sequentially test each of the M 
possible contexts. The selection of a proper schedule is 
based on prior knowledge about the amount of time 
required by the peak delay estimator to correctly locate a 
peak in different contexts. The first trial period tų should 
be related to the number of samples required to properly 
detect the global peak in an optimal context. The final trial 
period Ty reflects the amount of time willing to be spent 
operating in a worst-case context. The remaining 
trial periods can be uniformly spaced between these 
two limits. Intuitively, a schedule of non-decreasing 
trial periods ensures that the likelihood of selecting 
(yet-unattempted) contexts that can successfully lead to 
a correct estimate of the peak delay increases. At the same 
time, beginning the schedule with shorter waiting times 
escapes any suboptimal contexts earlier in the peak delay 
estimation process. Further analytical motivations in 
terms of the mean time to correctly locate the bulk delay 
are provided in Appendix B. 

In our implementation of this scheme, two counters 
Cinc and Cadec are used to maintain the number of samples 
that the peak tendency has either been increasing or 
decreasing during the current trial period Ty. If Cgec > Tp 
and jitter occurs in the peak delay estimate, then the filter 
and counters are reset to zero and a new context is 
attempted for the next trial period T,,;. If instead Cinc > Tx, 
i.e. peak delay is correctly estimated, then the counters are 
reset and the trial period is set back to its initial value 71. If 
all contexts have been attempted unsuccessfully, the 
process is repeated beginning with T>. 


3.4. Improved tracking (IT) 


Abrupt changes in an echo-path impulse response can 
be associated with a change in the bulk delay, or a phase 
roll, whereby the impulse response coefficients change 
signs. The tracking of a dispersive region after an abrupt 


change in the echo path can be seen as a competition 
between filter coefficient magnitudes. Indeed, when the 
partial Haar adaptive filter is initialized to zero, a peak’s 
magnitude only has to compete with the low-magnitude 
coefficient noise of the surrounding taps, which makes its 
detection easier. However, in the event of an abrupt 
change in the echo-path impulse response (as can occur 
on IP networks for example), a new peak might have to 
compete with the decreasing magnitude of an old peak. 
This can cause the peak delay estimator to take much 
longer to find the new peak. 

The proposed IT approach takes advantage of the 
following observation: After an abrupt change in the 
echo-path impulse response occurs, one finds that 
although the location of the new peak as well as its 
steady-state magnitude is unknown, the new steady-state 
magnitude of an old peak is approximately zero. Thus, if a 
significant decrease in the current peak’s magnitude is 
detected (which usually signals a change in the echo-path 
impulse response), the entire partial Haar filter is reset to 
v(n) = 0. This is a feasible solution because this filter is 
not directly being used to cancel echo. This allows the new 
peak to solely compete with the low-magnitude coeffi- 
cient noise of its neighbouring taps instead of the 
decreasing magnitude of the previous peak. As a result, 
the performance gains obtained in [14] for the stationary 
case can be extended to cases where abrupt changes in the 
echo-path impulse response occur. 

A reset is only deemed necessary when a decrease 
in magnitude is detected for a peak whose tendency 
has been in an increasing state during more than 
Tinc (not necessarily consecutive) iterations. To prevent 
the effects of false resets, the location of the peak before 
resetting the partial Haar adaptive filter, is stored and fed 
to the bulk delay unit. Only once the new peak tendency 
has been categorized as increasing for a total of Tinc 
samples, is the location of the new peak used to center the 
short time-domain filter. This way, if a reset is in fact 
necessary, the bulk delay unit will offset the short time- 
domain filter to the new peak location once it has been in 
an increasing state for a total of Tinc samples. If, however, a 
reset is unnecessary, then the old peak location is still 
used to center the short time-domain filter and the partial 
Haar adaptive filter will readapt to the old solution 
without affecting the time-domain adaptation. To avoid 
resetting the partial Haar adaptive filter too often, no new 
reset operation is attempted for at least Tp; iterations 
following a previous reset. In practice, Trs and Tinc can be 
selected as Trs ~ 0.25Tj,- © T1. The above IT approach can 
be easily integrated with the SCE approach in Section 3.3. 


3.5. Multiple echo (ME) cancellation 


The solution proposed here for ME cancellation with 
the PHDAF is based on a coordinator-multi-agent archi- 
tecture similar to [30]. The following assumptions are 
made: 


e The length of the dispersive regions is bounded by a 
known integer L in the time domain. 
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Fig. 6. Multiple echo cancellation system architecture. 


e The delay between successive echoes (or reflections) is 
larger than L, i.e. dispersive regions do not overlap. 

e The number of dispersive regions is upper bounded by 
a small integer np. 


These standard assumptions ensure that some advantages 
can be gained by using an adaptive solution tailored to the 
ME case; otherwise, the use of a full time-domain 
adaptive filter may be more appropriate. For instance, 
most network hybrid impulse responses range from 3 to 
12 ms, or a corresponding length of 24-96 samples at a 
sampling rate of 8kHz [20], which is typically much 
smaller than the maximum possible bulk delay. 

A block diagram of the proposed multiple echo-PHDAF 
(ME-PHDAF) is shown in Fig. 6. The partial Haar adaptive 
filter is partitioned into N/L sections, each of length L/M 
coefficients. Assigned to each region is an agent A; that 
runs a distributed form of the SCE algorithm that departs 
slightly from the case of a single dispersive region in 
Section 3.3. The length of each region is BL/M, obtained by 
extending the corresponding section of the partial Haar 
filter so that it overlaps with its left neighbour; the 
overlapping factor is typically set to f = 1.5. The PDM of 
agent i at time n is now given by 





Min(C; max; ĉ) 


PDM; = 1 — (24) 


Ci max 
Each agent sends its Cimin to the coordinator (see Fig. 6) 
where ĉis calculated as č = maxX; Cimin. The use of the min 
operator in (24) is necessary because ¢>C;max in some 
cases, Making the original PDM negative. This form of 
cooperating agents preserves the global characteristic of 
the PDM. 

If an agent requires a context switch, then only its 
corresponding section’s (not region) context should be 
changed. The polyphase decomposition of the RPHT from 
Section 3.1 comes in handy here, since each set of 
transformed input samples spanning a given section is 
readily available. The use of overlapping regions help to 
ensure that contiguous agents sharing a peak near their 
boundary have the same operating context, so that the 
orthogonality of the corresponding Haar basis vectors is 
preserved. Each agent is classified as active (a; = 1) or 
inactive (a; = 0), depending on its peak’s behavioural 


model. An agent becomes active if its peak tendency has 
been in an increasing state for more than Tin: samples; 
otherwise, it is deactivated. A maximum of np agents can 
be active at any given time, based on a first-in first-out 
mechanism. 

With information about the state of an agent, the 
central coordinator does one of two things: If one or more 
agents are active, the coordinator uses the peak locations 
of those agents to update only certain tap-weights of the 
time-domain adaptive filter centered on these locations. 
In effect, the input signal u(n) in Fig. 1 is time-shifted and 
applied as input to short adaptive filter sections of length 
L, represented by the shaded areas in Fig. 6. If none of the 
agents are active (i.e. during initialization of the algo- 
rithm), the coordinator uses a single-echo version of the 
SCE scheme to locate a single dispersive region and 
proceeds to update the coefficients of a time-domain filter 
of length L centered on this delay. 


3.6. Computational complexity 


For convenience, we shall refer to the PHDAF algorithm 
incorporating the proposed SCE and IT schemes of 
Sections 3.3 and 3.4 as the improved PHDAF (I-PHDAF), 
and to the PHDAF incorporating the multi-agent structure 
in Section 3.5 as the ME-PHDAF. We first note that the 
partial Haar transform used in these algorithms is 
computationally efficient as it does not require explicit 
multiplication operations. Indeed, for a given scale index j, 
the normalization factor &j is constant; it can be factored 
out of the matrix Hy and absorbed as part of subsequent 
operations on the transformed output vector. 

The bulk of the SCE and IT schemes’ computational 
complexity lies with the PTE. Table 3 shows the number of 
arithmetic operations per iteration required by the PTE, 
the SCE and the IT schemes, along with the required 
number of operations for an NLMS-based realization of 
the basic PHDAF structure in Fig. 1. We note that the 
arithmetic operations associated with the PTE are drasti- 
cally reduced since many bbm terms, such as v,;,(S A L) and 
i 4 UL), are zero throughout the PTS’s operation [21]. 
The logarithmic operations needed in (23) can also be 
eliminated since only the relative entropy is required. An 
alternate measure that preserves the desired relationship 
and requires only a single comparison operation is given 
by Hse = min(P” {S}, P“ {L}. Therefore, as an example, 
letting N = 1024, q = 256, and L = 128, the percentage 
increase in complexity when using the I-PHDAF algorithm, 
versus the basic PHDAF algorithm, is only 784, = 3.94%. 


Table 3 
Number of arithmetic operations per iteration—single dispersive region 
(I-PHDAF). 


PTE SCE IT PHDAF 
Add. 17 1 0 2q+2L+N/q+1 
Mult. 18 0 0 2q+2L+6 
Div. 3 0 0 1 
Comp. 11 5 6 0 
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Table 4 
Number of arithmetic operations per iteration—multiple dispersive 
regions (ME-PHDAF). 


Multi-agent PHDAF 
Add. 18N/L 2q+2nDL+N/q+1 
Mult. 18N/L 2q + 2nDL +6 
Div. 3N/L 1 
Comp. 18N/L —1 0 


The computational complexity of the ME-PHDAF is 
summarized in Table 4. Note that the complexity of the 
PHDAF component now includes provisions for the use of 
np short time-domain adaptive filters. The complexity of 
the multi-agent scheme is on the order of N/L times that 
of the basic PTE/SCE scheme in Table 3. To reduce the 
computational load, the ME-PHDAF can be programmed 
to allow only one agent to estimate its peak tendency 
every input sample, i.e. agents time-share the PTE. This 
also requires scaling the trial periods in the schedule t and 
relevant counters by N/L. 


4. Computer simulations 
4.1. Methodology 


The set of hybrid impulse responses used in the 
following simulations are taken from Annex D of the 
ITU-T G.168 Recommendation for digital network echo 
cancellers [20]. There are eight impulse responses m;,(n) 
(i= 1,2,...,8) with lengths L; that range from 64 to 128 at 
a Sampling rate of 8 kHz; it is assumed that m;(n) = O for 
n<0 and for n>L;. These impulse responses are time 
shifted as needed and scaled to obtain an echo return loss 
of 15 dB. The input samples u(n) are taken from a white 
zero-mean Gaussian process with unit variance. The error 
signal d(n) is obtained by filtering u(n) with the selected 
impulse response and adding a white Gaussian measure- 
ment noise v(n), uncorrelated with u(n) and of zero-mean 
and variance o2 = 107°3/10, The SNR is set to 30dB, 
unless stated otherwise. 

The length of the input data vector (i.e. maximum 
length of unknown echo path) is set to N= 1024. The 
PHDAF and proposed I-PHDAF and ME-PHDAF algorithms 
use a partial Haar transform of size q = 256 and the length 
of the dispersive region(s), i.e. short time-domain filter(s), 
is assumed to be L = 128. The partial Haar and short time- 
domain adaptive filters utilize the NLMS algorithm with a 
step-size u = 1.” In the proposed PHDAF algorithms, the 
initial context is set to 1, which corresponds to no shifting 
of the partial Haar basis vectors. A conventional NLMS 
algorithm (i.e. length N) is also used as a reference for 
comparison. 





> In the partial Haar NLMS update equation (4), ||u(n)||? is used 
instead of ||z(n)||? for the normalization step; this was found to make the 
algorithm more robust to the rank deficiency of the transform matrix Hy. 
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Fig. 7. Learning curves for ITU-T hybrid response ms(n) using a: (a) best 
and (b) worst-case bulk delay for the initial context used (SNR = 30 dB). 


4.2. Context escaping 


We use the following schedule of trial periods to test 
the different contexts: t = {150, 250,300,400}. Fig. 7(a) 
and (b) show the learning curves corresponding to an 
echo-path impulse response using ITU-T G.168 hybrid 
model m5(n — l) under the best and worst bulk delays 
(with respect to the initial context), respectively. The 
curves represent an ensemble of 200 runs for each 
simulation. 

For the best-case bulk delay (Fig. 7(a)), both the I- 
PHDAF and PHDAF show identical learning curves reach- 
ing steady state at around k = 750 compared to the NLMS 
which converges at around k = 5000. This is due to the 
fact that the NLMS adapts a far larger number of 
coefficients (1024 compared to 128). In fact, this result 
(i.e. when the best bulk delay is used for the initial 
context) may be viewed as a best-case performance for 
the proposed I-PHDAF (and also PHDAF) since it leads to 
the fastest possible detection of the global peak. Even if all 
the contexts are adapted in parallel, this would not lead to 
a faster convergence at the considered SNR level of 30 dB. 
For a worst-case bulk delay (Fig. 7(b)), the PHDAF never 
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Table 5 
Comparison of mean times and standard deviations to correctly estimate 
the peak delay for different SNRs. 


SNR (dB) I-PHDAF PHDAF 

Mean Std. Mean Std. 
30 91.5 75.4 121.1 203.0 
20 107.7 86.4 214.5 664.0 
15 167.4 138.3 362.7 1067.7 


10 421.4 387.1 531.7 1177.2 


seems to reach steady-state, while the proposed I-PHDAF 
converges much faster, nearly as well as in the optimal 
case, requiring about k = 1000 samples to converge. In 
effect, the I-PHDAF adds flexibility to the basic PHDAF and 
prevents it from “getting trapped” in suboptimal contexts. 
At the same time, the I-PHDAF is fairly stable, i.e. can 
remain “locked” on an optimal context. This behaviour 
may significantly increase the convergence speed of the 
echo canceller, as observed. 

In addition to learning curves, the mean time for each 
echo canceller to correctly estimate the location of a 
dispersive region was compared for different SNRs. Each 
row in Table 5 consists of the average and standard 
deviation over 500 independent runs. For each run, one of 
the eight ITU-T G.168 hybrid impulse responses is selected 
with equal probability, and a random bulk delay uni- 
formly distributed in the interval [0, 895] is included. In all 
cases, the proposed I-PHDAF algorithm finds the disper- 
sive region faster (the mean time is smaller) and more 
consistently (the standard deviation is much smaller). 
Although the PHDAF and the I-PHDAF display comparable 
mean times to convergence at very low SNR, the standard 
deviation of the PHDAF is three times larger. This reveals 
the robustness of using a fixed schedule t together with a 
PTE for different values of SNR. Of course, if the SNR does 
not change much over a specific channel, then schedules 
can be constructed specifically for those cases. 


4.3. Improved tracking 


For reasons of brevity, only one generic scenario will be 
considered here in which the bulk delay of the true echo- 
path impulse response abruptly changes from an optimal 
to a suboptimal delay with respect to the initial context. 
Specifically, the change is from m(n — 640) to ms5(n — 
322) at time k = 5200. The parameters of the proposed IT 
approach are chosen as Trs = 32 and Tinc = 128. 

The corresponding learning curves, averaged over 200 
independent runs are shown in Fig. 8. It can be seen that, 
initially, both echo cancellers converge optimally. How- 
ever, after the abrupt change in bulk delay, the I-PHDAF 
requires much fewer samples to converge. The learning 
curve of the PHDAF features a plateau-region extending 
for almost 2000 samples after the change. This corre- 
sponds to the amount of time it takes the new peak’s 
magnitude (which is very small since the new bulk delay 
is suboptimal) to exceed the decreasing magnitude of the 
old peak. 
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Fig. 8. Tracking behaviour of the I-PHDAF compared to the PHDAF and 
NLMS algorithm with a _ best-to-worst-case change in bulk delay 
(SNR = 30 dB). 


4.4. ME cancellation 


The ME-PHDAF algorithm of Section 3.5 is tested with 
N/L = 8 agents operating in parallel to classify and locate 
dispersive regions. The number of possible dispersive 
regions is set to np =3 and the parameter Tinc = 64 
samples. 

The learning curves in Fig. 9(a) and (b) first compare 
the ME-PHDAF’s performance to that of the PHDAF for the 
two cases of single echo path in Fig. 7. Each learning curve 
was averaged over 200 runs. For the best-case bulk delay 
in Fig. 9(a), the PHDAF is slightly faster than the ME- 
PHDAF, although both echo canceller’s require only 
around 1000 samples to converge. However, the distrib- 
uted form of the SCE algorithm allows the ME-PHDAF to 
converge again much faster in a worst-case bulk delay in 
Fig. 9(b). 

Fig. 10 shows the performance of the ME canceller 
when the number of dispersive regions in the true echo- 
path impulse response is 2. The learning curves are 
averaged over 600 independent runs consisting of bulk 
delays and dispersive regions randomly selected accord- 
ing to 


h(n) = gım;(n — 41)(N) + 8,m,(n — Ap), (25) 


where h(n) is the echo-path impulse response, g, and g» 
are fixed gains, i and j are independent equiprobable 
selections from {1,2,...,8}, 4; =320+1;, 42 = 640+ 12 
and r;,r2 are independent equiprobable in {0,1,2,3}. The 
IdNLMS filter represents the ideal case where exact prior 
knowledge about the location of these two dispersive 
regions is available, and only the corresponding time- 
domain filter coefficient are adapted while the remaining 
coefficients are set to zero. This approach, used as a 
benchmark, provides the best possible performance in this 
case. From the results in Fig. 10, we note that for k<500, 
the ME-PHDAF’s learning curve is slightly larger than the 
NLMS (only one peak is initially detected), after which its 
convergence speed increases, reaching the IdNLMS curve’s 
steady-state MSE at around k = 2500. 
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multiple echo canceller and the PHDAF: (a) best-case bulk delay; (b) 
worst-case bulk delay (SNR = 30 dB). 
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Fig. 10. Learning curves over 600 runs using two random hybrid impulse 
responses and random bulk delays (SNR = 30 dB). 


5. Summary and conclusion 


We have proposed feasible solutions to improve the 
performance of the partial Haar dual adaptive filter (PHDAF) 
in practical applications. These includes: (1) alleviating the 
dependence of the PHDAF’s performance on the bulk delay 
of the echo-path impulse response; (2) improving the 
tracking performance of the PHDAF in response to abrupt 
changes in the echo path; and (3) extending the original 
PHDAF structure to support the cancellation of MEs. 

The proposed algorithmic solutions exploit the poly- 
phase representation of the RPHT and make use of a peak 
tendency estimator (PTE) based on Dezert-Smarandache 
theory (DSmT) and fuzzy inference. The PTE monitors the 
time evolution of the Haar-domain peak coefficient’s 
discernibility, and categorizes it as being in an increasing 
or a decreasing state. In the proposed IPHDAF algorithm, 
this information is used to decide if the current operating 
context (i.e. polyphase component of the Haar transform) 
should be changed, as well as to detect if an abrupt change 
in the echo-path impulse response has occurred. To deal 
with the case of multiple dispersive regions, an ME-PHDAF 
algorithm is also proposed which uses a distributed (i.e. 
multi-agent) form of the above PTE-based approach. 

The IPHDAF and ME-PHDAF algorithms were evaluated 
in terms of their mean-square error (MSE) curves as well 
as their mean time to properly locate dispersive regions 
for different SNR and channel conditions. Results show 
that a significantly enhanced performance can be ob- 
tained using the proposed algorithms at a minimal 
increase in computational cost when compared to the 
original PHDAF algorithm. Finally, we note that the 
techniques proposed in this paper for the PHDAF are not 
conceptually limited to the Haar transform and could be 
applied to other types of wavelets as well. 


Appendix A. Theoretical motivation behind the PDM 


In this appendix, we further motivate the choice of the 
PDM (13) on theoretical grounds. We begin by deriving an 
expression for Pc, defined as the probability of correctly 
locating the peak magnitude of the partial Haar adaptive 
filter coefficient vector v, where the dependence on 
discrete-time n is omitted for convenience. To this end, 
we generalize the development in [14] to the case of an 
arbitrary mean vector. Specifically, let v = [v1,...,vg] and 
assume that the individual weights v; are statistically 
independent random variables with mean c; and prob- 
ability density functions (pdf) f,(u) = (u — ci), derived 
from a common symmetric pdf (u) with zero-mean, i.e. 
denote the corresponding index; without loss in general- 
ity, assume that Ci, = Cmax (positive peak). Then, Pc can be 
expressed as 


Pc = l Pr{|vj|<|vj,| for all i#io|Vi = u}p(u — Cmax) du 
(A.1) 


= CAM (uj; — ci) cule —Cmax)du. (A.2) 
—& (izio Y lil< lul 
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Next, we derive a lower bound and an upper bound for Pc. 
Let .4,,...,-%~ define a partition of the set J = {ice N: 
1<i<q andi¥ig} and let p, denotes the number of 
elements in p. For each subset .%;, define b= 
Minjes, (lci) and dy = MaXjez, (|Cj|). Due to the symmetry 
of ow), it can be verified that PL <Pc <P where 


oT Pk 
= ~~ — 
Pe = [ I] ff. <|ul P(Uk dy) du} plu Cmax) du 


CO k=1 


(A.3) 


and P¢ is given by a similar expression with the numbers 
d replaced by b;, respectively. We note that this bound on 
Pc becomes tighter as the size p, of the subsets p is 
reduced. 

Define Cmin = MINg_1.... x ([dx|), i.e. the minimum among 
the maximum peaks of each of the intervals p. We 
assume that each subinterval x, contains at least one 
coefficient with peak magnitude less than Cmin, Le. 
bk <Cmin, Which is reasonable for larger subintervals. Then 
we also have P- <P <P" where 


ore) q-1 
c = J if P(V — Cin) dv) P(U — Cmax)du (A.4) 
—oo (J |v|<|ul 
now provides an approximation to Pc. It can be shown 
that: P- is an increasing function of Cmax; Pc is a 
decreasing function of Cmn; P71 in the limit 
Cmax >Cmin} and P- ~ O in the limit Cmin = Cmax and q>1. 
This behaviour is captured by the ideal PDM 


PDM? = 1 — Cmin/Cmax- (A.5) 


In this respect, Eq. (13) represents a practical attempt at 
estimating the unknown quantity PDM? based on the 
available data, i.e. the filter coefficients v;(n) at time n, and 
the use of only three subintervals. 


Appendix B. Motivation behind using a schedule of trial 
periods 


The motivation for using a schedule of trial periods is 
embedded in system reliability theory [31]. We assume 
the existence of M systems Y4,...,YmĮm that provide the 
solution to the same problem. The time for each system to 
produce the solution is a random variable governed by an 
exponential cumulative distribution F(t; 2) = 1 — e~“' for 
i=1,...,M, with mean 4; = 1/2; and variance o? = 1/27, 
where 4; is the rate parameter for system i. A randomly 
ordered selection of the M systems is made, as repre- 
sented by a 1-to-1 permutation z of the integers 
{1,...,M}. Then, beginning with system Y,;, up to Yzy, 
each system is tested sequentially for a period of time 
T1,...,T7m, respectively. Given knowledge of the different 
rates j;, and assuming that the systems operate indepen- 
dently, the problem of interest here is to choose the trial 
times T;, or equivalently the system switching times 
Si =71+---+17;, so that the total time spent is searching 
a solution is minimized on a Statistical basis under the 
constraint of a fixed total trial time T = sy. The problem is 
made difficult by the fact that the user does not know the 
permutation 7, i.e. which system is being used in a given 
trial interval. 


Let random variable X denote the total time needed to 
solve the problem. The conditional CDF of X given a 
permutation z is shown to be 


Fx(t|70) =l- ap(T)e “7K ES), Sk-1 S C<Sx, (B.1) 


where we define (mT) = k-i Ant for k>1, um) = 1 


and so = 0. For t<0, we have Fx(t|z) = 0 while for t>sSm, 
we assume Fy(t|z) = 1, i.e. the probability of not obtaining 
a solution within a total time sy is negligible. The 
conditional expectation of X can be obtained as 


M 


E[X|z] = [a — Fx(t\7)) dt = SH 


k=1 “Tk 


—e 4%), (B.2) 


Finally, assuming that the different possible permutations 
are equiprobable, the expected total solution time is 
obtained as 


1 
E[X] = vie EXI], (B.3) 


where the summation is over all possible permutations z. 

While an exhaustive analysis of the above expression 
for E[X] is beyond the scope of the paper, we can 
nevertheless obtain interesting insights by considering 
M = 2, in which case (B.3) takes the simplified form 





me 1 1 1 —A, T1 A2 —A2T2 
1 1 —A2T Ay —A1T2 


In the case of equal rates, i.e. 2; = A> = å, this expression 
simplifies to E[X] = (1/2)(1 — e~*") which is independent 
of the switching time tT; = s1. However, in the case 1; 4/2, 
the situation is different, as revealed by the plots of E[X] 
versus Tı in Fig. B1, obtained under the constraint t; + 
T2 = T for different values of the rate parameters 2, and 
A2. The value of tT; minimizing E[X] in these plots can be 
obtained by applying the method of Lagrange multipliers 
and satisfies (0<T, <T): 


Aye-*2(1 — eT BOD) = Agen — eT), (B.5) 
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Fig. B1. Plot of E[X]/T vs. t; /T for three pairs of (21T, 22T). 
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We generally find that for a wide range of parameter 
values under consideration, the optimum switching time 
occurs at tT; ~ 0.3T. Clearly, this approach can be general- 
ized to the case of M>2 intervals. 
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